home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Language/OS - Multiplatform Resource Library
/
LANGUAGE OS.iso
/
cpp_libs
/
vector.lha
/
vector
/
cmatrix.h
< prev
next >
Wrap
C/C++ Source or Header
|
1991-11-23
|
6KB
|
244 lines
#ifndef _CMATRIX_H
#define _CMATRIX_H
#include "cvector.h"
#include "enum.h"
// Nest axis enumeration so it must be qualified by e.g. Axis::x
class Axis {
public:
enum axis { x, y, z };
};
#define CartesianMatrixDeclare(T) \
\
typedef T CMatType(T)[3][4]; \
\
class CMat(T) { \
protected: \
T m[3][4]; \
\
public: \
const T *operator()(int i) const { \
return m[i]; \
} \
T *operator[](int i) { return m[i]; } \
operator CMatType(T)() { return m; } \
\
CMat(T) &zero(); \
CMat(T) &identity(); \
}; \
\
/* V = M * V (column vectors assumed) */ \
CVec(T) operator*(const CMat(T) &m, const CVec(T) &v); \
\
/* M = M * M */ \
CMat(T) operator*(const CMat(T) &a, const CMat(T) &b); \
\
/* Geometric transformations */ \
CMat(T) rotation_matrix(Enum(Axis,axis) a, T alpha); \
CMat(T) rotation_matrix(const CVec(T) &axis_, T alpha); \
CMat(T) translation_matrix(const CVec(T) &t); \
CMat(T) scale_matrix(T s); \
CMat(T) scale_matrix(const CVec(T) &s); \
\
ostream &operator<<(ostream &o, const CMat(T) &m);
#define CartesianMatrixImplement(T) \
\
static const int \
X = 0, \
Y = 1, \
Z = 2, \
W = 3; \
\
CMat(T) &CMat(T)::zero() { \
for (int i = X; i <= Z; i++) \
for (int j = X; j <= W; j++) \
m[i][j] = 0; \
return *this; \
} \
\
CMat(T) &CMat(T)::identity() { \
for (int i = X; i <= Z; i++) \
for (int j = X; j <= W; j++) \
m[i][j] = (i == j); \
return *this; \
} \
\
CVec(T) operator*(const CMat(T) &m, const CVec(T) &v) { \
return Vector( \
m(X)[X]*v(X)+m(X)[Y]*v(Y)+m(X)[Z]*v(Z)+m(X)[W], \
m(Y)[X]*v(X)+m(Y)[Y]*v(Y)+m(Y)[Z]*v(Z)+m(Y)[W], \
m(Z)[X]*v(X)+m(Z)[Y]*v(Y)+m(Z)[Z]*v(Z)+m(Z)[W]);\
} \
\
CMat(T) operator*(const CMat(T) &a, const CMat(T) &b) { \
CMat(T) res; \
\
for (int i = X ; i <= Z; i++) { \
res[i][X] = a(i)[X] * b(X)[X] + \
a(i)[Y] * b(Y)[X] + \
a(i)[Z] * b(Z)[X]; \
res[i][Y] = a(i)[X] * b(X)[Y] + \
a(i)[Y] * b(Y)[Y] + \
a(i)[Z] * b(Z)[Y]; \
res[i][Z] = a(i)[X] * b(X)[Z] + \
a(i)[Y] * b(Y)[Z] + \
a(i)[Z] * b(Z)[Z]; \
res[i][W] = a(i)[X] * b(X)[W] + \
a(i)[Y] * b(Y)[W] + \
a(i)[Z] * b(Z)[W] + a(i)[W]; \
} \
\
return res; \
} \
\
static void sincos(T &cval, T &sval, T alpha) { \
/* Tolerance before assuming the rotation \
* is aligned with axes \
*/ \
const T tolerance = 1e-10; \
\
cval = cos(alpha); \
sval = sin(alpha); \
\
if (cval > 1 - tolerance) { \
cval = 1; \
sval = 0; \
} else if (cval < -1 + tolerance) { \
cval = -1; \
sval = 0; \
} \
\
if (sval > 1 - tolerance) { \
cval = 0; \
sval = 1; \
} else if (sval < -1 + tolerance) { \
cval = 0; \
sval = -1; \
} \
} \
\
/* Create a rotation matrix of alpha radians around \
* specified axis (x,y,z). \
*/ \
CMat(T) rotation_matrix(Enum(Axis,axis) a, T alpha) { \
T ca, sa; \
sincos(ca, sa, alpha); \
CMat(T) r; \
\
r.zero(); \
r[W][W] = 1; \
cout << "rotation_matrix: x = " \
<< Axis::x << ' ' << (int)Axis::x \
<< " a = " << a << ' ' << (int)a << endl; \
\
switch (a) { \
case Axis::x: \
r[X][X] = 1; \
r[Y][Y] = ca; r[Y][Z] =-sa; \
r[Z][Y] = sa; r[Z][Z] = ca; \
\
break; \
case Axis::y: \
r[X][X] = ca; r[X][Z] = sa; \
r[Y][Y] = 1; \
r[Z][X] =-sa; r[Z][Z] = ca; \
\
break; \
case Axis::z: \
r[X][X] = ca; r[X][Y] =-sa; \
r[Y][X] = sa; r[Y][Y] = ca; \
r[Z][Z] = 1; \
\
break; \
default: \
cerr << "rotation_matrix: unknown axis " \
<< (int)a << ", returning identity" \
<< endl; \
r.identity(); \
break; \
} \
\
return r; \
} \
\
/* Create a rotation matrix around specified vector */ \
CMat(T) rotation_matrix(const CVec(T) &axis_, T alpha) {\
T ca, sa; \
sincos(ca, sa, alpha); \
CMat(T) c, s, r; \
\
/* rotation axis must be normalized */ \
CVec(T) a = axis_; \
a.normalize(); \
\
/* matrix effecting P' = (1 - cos alpha) (P dot A) A\
* Rij = (1 - cos alpha) Ai Aj \
*/ \
for (int i = X; i <= Z; i++) \
for (int j = X; j <= Z; j++) \
c[i][j] = (1 - ca) * a(i) * a(j); \
\
/* matrix effecting P' = (cos alpha) P + (sin alpha) P ^ A */ \
s.zero(); \
s[X][X] = ca; s[X][Y] =-sa * a(Z); s[X][Z] = sa * a(Y); \
s[Y][X] = sa * a(Z); s[Y][Y] = ca; s[Y][Z] =-sa * a(X); \
s[Z][X] =-sa * a(Y); s[Z][Y] = sa * a(X); s[Z][Z] = ca; \
\
/* overall rotation is sum of s and c */ \
for (i = X; i <= Z; i++) \
for (j = X; j <= Z; j++) \
r[i][j] = s[i][j] + c[i][j]; \
r[X][W] = r[Y][W] = r[Z][W] = 0; \
\
return r; \
} \
\
/* Create a translation matrix by vector t */ \
CMat(T) translation_matrix(const CVec(T) &t) { \
CMat(T) r; \
\
r.identity(); \
r[X][W] = t(X); \
r[Y][W] = t(Y); \
r[Z][W] = t(Z); \
\
return r; \
} \
\
/* Create an isotropic scaling matrix */ \
CMat(T) scale_matrix(T s) { \
return scale_matrix(CVec(T)(s,s,s)); \
} \
\
/* Create an anisotropic scaling matrix */ \
CMat(T) scale_matrix(const CVec(T) &s) { \
CMat(T) r; \
\
r.zero(); \
r[X][X] = s(X); \
r[Y][Y] = s(Y); \
r[Z][Z] = s(Z); \
\
return r; \
} \
\
ostream &operator<<(ostream &s, const CMat(T) &m) { \
for (int i = X; i <= Z; i++) \
s << "\t" << "( " \
<< m(i)[X] << " " \
<< m(i)[Y] << " " \
<< m(i)[Z] << " " \
<< m(i)[W] << " )\n"; \
\
return s; \
}
CartesianMatrixDeclare(float);
typedef CMat(float) Matrix;
#endif /*_CMATRIX_H*/